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Galaxy stability within a self-interacting dark matter halo 
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ABSTRACT 

This paper investigates spheroidal galaxies comprising a self-interacting dark matter 
halo (SIDM) plus de Vaucoulcurs stellar distribution. These are coupled only via their 
shared gravitational field, which is computed consistently from the density profiles. 
Assuming conservation of mass, momentum and angular momentum, perturbation 
analyses reveal the galaxy's response to radial disturbance. The modes depend on 
fundamental dark matter properties, the stellar mass, and the halo's mass and radius. 
The coupling of stars and dark matter stabilises some haloes that would be unstable 
as one-fluid models. However the centrally densest haloes are unstable, causing radial 
flows of SIDM and stars (sometimes in opposite directions). Depending on the dark 
microphysics, some highly diffuse haloes are also unstable. Unstable galaxies might 
shed their outskirts or collapse. Observed elliptical galaxies appear to exist in the safe 
domain. Halo pulsations are possible. The innermost node of SIDM waves may occur 
within ten half-light radii. Induced stellar ripples may also occur at detectable radii 
if higher overtones are excited. If any SIDM exists, observational skotoseismology of 
galaxies could probe DM physics, measure the sizes of specific systems, and perhaps 
help explain peculiar objects (e.g. some shell galaxies, and the growth of red nuggets). 



Key words: dark matter — galaxies: elliptical and lenticular, cD — galaxies: kine- 
matics and dynamics — galaxies: structure 



1 INTRODUCTION 

For decades it has been observed that visible matter in 
galaxies and clusters is insufficient to bind tliem gravita- 
tionally (e.g. Zwicky 1937; Freeman 1970; Roberts & Rots 
1973). The extra gravity of invisible "dark matter" (DM) is 
invoked, but its properties are uncertain, with many parti- 
cle physics conjectures offering candidates. On Earth, DM 
particle detection experiments differ in their findings: some 
report annually modulated signals (e.g. Bernabei et al. 2008; 
Aalseth et al. 2011a,b; Angloher et al. 2012); others see 
nothing above their expected background (e.g. Angle et al. 
2008; Aprile et al. 2010; Akerib et al. 2006; Ahmed et al. 
2009). Reconciling them may require extra physics: inelastic 
or composite forms of dark matter, or extra internal degrees 
of freedom (e.g. Smith et al. 2001; Chang et al. 2009; Alves 
et al. 2010; Kaplan et al. 2010, 2011) 

Large-scale cosmic filaments and voids emerge in A^- 
body computer simulations using gravity alone (e.g. Frenk, 
White & Davis 1983; Melott et al. 1983; Springel, Frenk & 
White 2006) but there is observational tension concerning 
smaller, self-bound systems. The predicted abundances of 
dwarf galaxies exceed those seen in many types of environ- 
ments (e.g. Klypin et al. 1999; Moore et al. 1999; D'Onghia 
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& Lake 2004; Tikhonov & Klypin 2009; Carlberg, Sullivan 
& Le Borgne 2009; Zavala et al. 2009; Peebles & Nusser 
2010). The biggest simulated subhaloes appear too massive 
to match real satellite galaxies (Boylan-Kolchin, Bullock & 
Kaplinghat 2011, 2012). CoUisionless DM haloes develop sin- 
gular central density cusps (e.g. Gurevich & Zybin 1988; Du- 
binski & Carlberg 1991; Navarro, Frenk & White 1996). On 
the contrary, observations of diverse galaxy classes indicate 
or allow a dark density profile that is smooth and nearly 
uniform within a kpc-scale inner core, (e.g. Flores & Pri- 
mack 1994; Moore 1994; Burkert 1995; Salucci & Burkert 
2000; Kelson et al. 2002; Kleyna et al. 2003; Goerdt et al. 
2006; Gentile et al. 2004; de Blok 2005; Thomas et al. 2005; 
Kuzio de Naray et al. 2006; Gilmore et al. 2007; Weijmans 
et al. 2008; Oh et al. 2008; Nagino & Matsushita 2009; In- 
oue 2009; de Blok 2010; Pu et al. 2010; Murphy, Gebhardt 
& Adams 2011; Memola, Salucci & Babic 2011; Walker & 
Peharrubia 2011; Jardel & Gebhardt 2012; Agnello & Evans 
2012; Salucci et al. 2012; Lora et al. 2012). Galaxy cluster 
haloes should be centrally compressed by the steep density 
profile of cooling gas, but nontheless many examples are def- 
initely or possibly cored (e.g. Sand, Treu & Ellis 2002; Sand 
et al. 2008; Ettori et al. 2002; Halkola, Seitz & PanneUa 
2006; Halkola et al. 2008; Voigt & Fabian 2006; Rzepecki 
et al. 2007; Zitrin & Broadhurst 2009; Newman et al. 2009, 
2011; Richtler et al. 2011). 
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To remedy (or defer) these DM crises, many astro- 
physicists invoke energetic "feedback" by starbursts or ac- 
tive galaxies: explosions and outflows violent enough to flat- 
ten the dark cusps and denude small subhaloes of baryons 
(e.g. Navarro, Eke & Frenk 1996; Gnedin & Zhao 2002; 
Mashchenko, Couchman & Wadsley 2006; Peirani, Kay & 
Silk 2008; Governato et al. 2010). These events would have 
happened early enough to be difhcult to observe today. The 
details would involve radiative shocks in magnetised multi- 
phcise interstellar media that are thermally and dynamically 
unstable. The relevant scales are numerically unresolvable 
in present simulations. Results depend on phenomenologi- 
cal recipes for sub-grid processes, and the implementations 
are more varied and elaborate than any DM physics. For rea- 
sons of resolution alone, the proposition that feedback cures 
the DM crises cannot be decisively evaluated for at least 
another decade. Regardless, feedback may be energetically 
incapable of removing cusps whilst leaving realistic stellar 
populations in dwarf galaxies (Pefiarrubia et al. 2012). 

An alternative approach extends standard particle 
physics phenomenologies: taking the dark matter crises as 
signatures of dynamically important physics in the dark 
sector. In "SIDM" theories, dark matter is self-interacting 
yet electromagnetically neutral. Pressure support or ther- 
mal conduction could prevent cusps, while drag and abla- 
tion reduce subhaloes (e.g. Flores & Primack 1994; Moore 
1994; Spergel & Steinhardt 2000; Firmani et al. 2000; Vogels- 
berger, Zavala & Loeb 2012; Rocha et al. 2012). SIDM inter- 
actions entail direct inter-particle scattering (with specific 
crossection ?) possibly enhanced by short-range Yukawa- 
type forces, or a long-range "dark electromagnetism" (e.g. 
Ackerman et al. 2009; Loeb & Weiner 2011). CoUisionless 
plasma supported by "dark magnetism" is conceivable. Al- 
ternatively, the SIDM might consist of boson condensates 
of a scalar field (e.g. Peebles 2000). Plausible early-universe 
evolution and cosmic structures emerge in cosmogonies with 
fluid-like SIDM (e.g. Moore et al. 2000; Stiele, Boeckel & 
Schaffner-Bielich 2010, 2011) or scalar field dark matter 
(e.g. Matos & Arturo Urena-Lopez 2001; Suarez & Matos 
2011). There were early fears that SIDM haloes would form 
sharper cusps via gravothermal collapse (e.g. Burkert 2000; 
Kochanek & White 2000). Later analyses put this beyond 
a Hubble time, and found two safe regimes that suppress 
thermal conductivity: weaker scattering (rarely per orbit); 
and stronger (fiuid-like) cases with tiny mean free paths 
(Balberg, Shapiro & Inagaki 2002; Ahn & Shapiro 2005; 
Koda & Shapiro 2011). We will refer to the <; < Icm^g"^ 
and c > 200cm^g~^ regimes as "weak-SIDM" and "strong- 
SIDM" respectively. 

The ubiquity of dark cores supports SIDM, but there 
have been circumstantial counter-arguments. It was argued 
that observed halo ellipticities imply low coUisionality (e.g 
Miralda-Escude 2002; Buote et al. 2002). Such limits can- 
not apply to systems that are elliptical due to tidal distur- 
bance, rotation or non-radial pulsations. Newer SIDM sim- 
ulations reject the spherical prediction and imply a looser 
limit c < Icm^g"^ (Peter et al. 2012). Other papers in- 
ferred tight limits on SIDM coUisionality, from the growth of 
supermassive black holes (? < 0.02cm^g~^: Ostriker 2000; 
Hennawi & Ostriker 2002) and the survival of cluster mem- 
bers (< 0.35 or > 10''cm2g-\ Gnedin & Ostriker 2001). 
However those investigations were in a way self-fulfilling, be- 



cause they assumed cuspy initial conditions, with a strong 
thermal inversion at the centre. This choice was inconsistent 
with cored profiles expected from SIDM. 

Gravitational lensing models of merging galaxy clus- 
ters give conflicting signs about DM physics: sometimes the 
inferred DM seems collisionless, following galaxies and sep- 
arating from gas (c < 1.25cm^g~^, e.g. Clowe et al. 2006; 
Randall et al. 2008). In other systems DM seems collisional, 
following the gas or separating from galaxies (e.g. Mahdavi 
et al. 2007; Jee et al. 2012; Williams & Saha 2011). The 
derived limits depend on assumed projected geometries and 
histories of unique systems; these are each open to specific 
counter-explanations and rebuttals. For instance if the ac- 
claimed "Bullet Cluster" were truly a direct-hit merger then 
ram pressure and shocks would quench star formation. This 
is not observed (Chung et al. 2009). 

Saxton & Wu (2008) modelled galaxy clusters in which 
the halo could possess extra heat capacity (e.g. composite 
SIDM particles). Stationarity in the presence of gas infiows 
requires a central accretor exceeding some minimum mass. 
For the range of DM properties giving plausible cluster and 
core profiles, a part of the inner halo teeters towards local 
gravitational instability, implying likely collapse to form su- 
permassive black holes of realistic sizes. Saxton & Ferreras 
(2010) modelled elliptical galaxies comprising a polytropic 
halo plus coterminous collisionless stars. Observed stellar 
kinematics of elliptical galaxies were fitted by the same class 
of DM that gives the most realistic clusters. 

This paper is a further theoretical study of spheroidal 
galaxies. Their stability properties will imply constraints on 
the end-points of galaxy evolution. If galaxies do possess 
SIDM haloes then their oscillatory modes will open the fu- 
ture possibility of skotoseismology — acoustic diagnosis of 
dark matter haloes that pulsate. The paper's organisation is 
as follows: §2 formulates the stationary model and perturba- 
tions; §3 describes the oscillatory solutions and instabilities; 
§4 discusses applications and interpretations; §5 concludes. 



2 FORMULATION 

This paper investigates spherical galaxies comprising self- 
interacting dark matter and collisionless stars in near- 
equilibrium configurations. The SIDM halo is polytropic, 
adiabatic, non-singular and radially finite (§2.1). The sys- 
tem is isolated: unconfined and unexcited by any external 
medium. Stellar orbits are isotropically distributed, in an 
observationally motivated density profile (§2.2). Oscillations 
and instabilities are calculated in response to radial per- 
turbations (§2.4-§2.6). Mass displacements are non-uniform, 
and at any location the two constituents may be displaced 
with differing amplitudes and phases. The perturbed vari- 
ables are non-singular at the origin, and the outer boundary 
of the halo is free (§2.7). 

2.1 dark matter halo profile 

For each mode in which ideal collisional gas particles are free 
to move, the energy per particle is |feT = ^fJ-o-^ where k is 
Boltzmann's constant, T is temperature, jj, is the particle 
mass and is a characteristic (one-dimensional) velocity 
dispersion. If SIDM gas has F effective thermal degrees of 
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freedom that are locally in equipartition (due to coUisional- 
ity or collective effects), then the specific heat capacities at 
constant volume and pressure are cv ~ \Fk and cp = cv+fc. 
The ratio of specific heats 

7=^ = 1+|. (1) 



F 



From the first law of thermodynamics, adiabatic processes 
imply a polytropic equation of state for the dark pressure: 

P,^p,al = spl, (2) 

or ( equivalent ly) the mass density is = Q o"d ^ where s is a 
pseudo-entropy, and Q = s~^''^ is a generalised phase-space 
density. In adiabatic processes, s and Q remain constant. 

The micro-physical meanings of deserve further re- 
marks. Point-like particles moving in three spatial dimen- 
sions have _F = 3, due to translational kinetic energy. Com- 
posite particle substructure raises F (e.g. rotational and vi- 
brational modes). A diatomic gas has _F = 5 because two 
parts are moving in 3D space, at a fixed separation. In- 
compressible fluids are more constrained and F = 0. For a 
relativistic or radiation dominated plasma, F = 6. If par- 
ticles were massive enough (short wavelength) to traverse 
hidden spatial dimensions, then extra translational modes 
would raise F. An isothermal medium has F oo. 

Polytropic fluids also emerge in some non-classical theo- 
ries of dark matter. Scalar fleld and boson condensate mod- 
els usually imply F = 2, but possibly other F values de- 
pending on an index in the interaction potential (e.g. Good- 
man 2000; Peebles 2000; Arbey, Lesgourgues & Salati 2003; 
Bohmer & Harko 2007; Harko 2011a,b; Chavanis & Delfini 
2011; Robles & Matos 2012). This paper is agnostic about 
which particle theory causes the value of F, which should 
ultimately be constrained by observations. 

The radial structure of a spherical halo is determined 
by its momentum equation, 

(3) 



d^r dP 
Pd — + 



di2 ■ dr 

The stellar mass enclosed at radius r is m.^ = m^^(< r) (see 
§2.2). The corresponding enclosed dark mass is obtained by 
numerical integration of 



(4) 



(5) 



dm. , 
dr ' 

In hydrostatic equilibrium, the condition (3) becomes 

d^' 2 ^(m^-fynj 

dr F + 2 r2 

Neglecting any compact central mass (e.g. black hole), 
a non-singular boundary condition applies at the origin 
(dcr^/dr|r=o ~ 0). The equilibrium dark halo is adiabatic, 
and its pseudo-entropy s is spatially uniform. This may oc- 
cur in nature if a galaxy halo is well mixed, or if s is set by 
universal constants of the SIDM particle. 

For —2 < _F < 10 a halo of finite mass truncates 
(/9j — > 0) at a finite outer radius R. For other F, any finite 
mass equilibrium model has infinite radius: it boils off into 
the void, or fails to condense from the cosmic background 
in the first place. The high limit {F — 10) is a generalised 
Plummer (1911) model, with infinite radius and finite mass. 
Haloes may also seem radially infinite in some SIDM simula- 
tions with _F = 3 (e.g. Rocha et al. 2012). This may be due to 



ongoing supersonic infall of neighbouring objects, or because 
the scattering is weak enough to leave the outskirts collision- 
less. If the SIDM "strength" is increased (greater scattering 
1,, or long-range "dark magnetism" effects) then polytropic 
behaviour affects more of the halo. Beyond some threshold 
strength, a definite surface must occur. Rocha et al. (2012) 
note a temperature drop at large radii in their largest c sim- 
ulation. This is attributable to incipient truncation. 

The value of F determines the natural ratio of core to 
system radii: for greater F the core is smaller (e.g. compare 
cases in Figure 1). If dark mass dominates, the circular orbit 
velocity profile {V = ^ Gm^/r) peaks at a radius i?o which 
occurs at Ro/R « 0.750, 0.0782 when F = 3 and F = 9 
respectively. In these cases, the dark density index drops 
to dln/9j/dlnr = -1 at radn R^/R ^ 0.379,0.0285; and 
index —2 a.t R2/R ~ 0.520, 0.0459. In other words, a strong- 
SIDM core with F = 9 is radially smaller than a tenth of 
a F = 3 core. Appendix A lists standard radii for other 
F cases. Adding stars distorts these ratios but the general 
trend with F persists. Figure 1 (right panel) shows how a 
stellar component can shrink and steepen the dark core. For 
comparison, the "scale radius" of simulated CDM haloes 
(Navarro, Frenk & White 1996) is R2, and Ro ^ 2.163i?2. 
For the Burkert (1995) cored profile with scale radius Vh, the 
key radn are Ri « 0.657rb, R2 ~ 1.521rb, and Ro « 3.244rb. 
Simulated SIDM cores fitted by Rocha et al. (2012) (rb « 
0.7 R2) appear consistent with Burkert's parameterisation. 
The ratios R\ : R2 : Ro are similar for Burkert models and 
F = 3 polytropes, but less similar for other F. 

These ideal core sizes are upper limits. In coUisional 
SIDM, weaker ? delays the formation of full-sized cores. 
Current SIDM simulations focus on the F = 3 case, which 
predicts cores larger than those observed, unless the scat- 
tering crossection is kept small (? < Icm'^g"^, e.g. Yoshida 
et al. 2000; Arabadjis, Bautz & Garmire 2002; Katgert, Bi- 
viano & Mazure 2004; Vogelsberger, Zavala & Loeb 2012; 
Rocha et al. 2012). The alternative is to vary F: the range 
of 7 < F < 10 results in galaxy clusters with realistic ~ 10^- 
lO^kpc cores (Saxton & Wu 2008). In principle, haloes with 
F > 7 could produce cores of realistic radius for arbitrarily 
large Such F values also provide Ro/R small enough to 
be consistent with observations of large flat regions in disc 
galaxies' rotation curves. The choice of F > 7 is therefore 
an astrophysically preferable subdomain, but for the sake of 
generality this paper considers a wider interval, 2 ^ F ^ 9. 



2.2 stellar mass profile 

The assumed stationary distribution of the stars resembles 
the empirical Sersic (1968) form, where the surface density 
at projected radius X is 



^■k{X) oc exp [-6(X/7?c)^''"] 



(6) 



and n is the "Sersic index." The "effective radius" Ro en- 
circles half the starlight. Assuming that the stellar mass-to- 
light ratio is constant, there exist exact deprojections of (6) 
that express the 3D density profile (e.g. Mazure & Capelato 
2002) , but these involve exotic special functions. We will in- 
stead use the approximate deprojection of Prugniel & Simien 
(1997), with density and mass profiles 
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radius r / 



Figure 1. Radial density profiles of equal-mass haloes with F = 2, . . . , 9 (coloured) truncating at the same radius {R = lOO-Ro). The left 
panel shows pure polytropes; the right panel includes a stellar profile (§2.2, "-A-") massing 1/24 times the halo. The halo density index is 
— 1 and —2 at the O and symbols (radii iJi and -R2) just outside the uniform core. Ticks locate the peak circular orbit velocity {Ra)- 



p-k = pc exp \ h-h (r/Rc) 



and 



(7) 



m*(< r) = 47rn6"(''-^'e' p.R^ T [n(3 - p), b (r/i?e)'/"] ,(§) 

where T{a, u) is the lower incomplete gamma function. The 
index n determines b and p (Lima Neto, Gerbal & Marquez 
1999; Ciotti & Bertin 1999; Marquez et al. 2000). Adopting 
n = 4 to match the classic de Vaucouleurs (1948) profile, 
the values are 6 « 7.669 and p ~ 0.8510. For this standard 
choice, the total stellar mass is M^f = 45.79pci??, and that 
within the half-light radius is m^{lRc) — IQ.OSpcRi- 

The frequency scale of the stars will be abbreviated as 

OJc = \/G Pc ■ (9) 

Given the properties of large elliptical galaxy like M87, 
Pe ~ 0.04 7710 pc and cjc ^ 0.01 Myr"^ (Murphy, Geb- 
hardt & Adams 2011). For a compact elliptical such as M32, 
Pc ~ 14m0pc~^ and ojc ~ 0.3 Myr"^ (Mateo 1998; Choi, 
Guhathakurta & Johnston 2002). For "red nuggets" at high 
redshift, pc > 3m0pc" ^ and ujc < 0.05 Myr"^ (e.g. Toft 
et al. 2012; van de Sande et al. 2011). If stars were the only 
component present, the keplerian period at r = IR^ would 
be ~ 1.440/a;e, and the angular frequency Hc — 4.363a;e. 
It is convenient to write the density gradient index as 



d In p^ 



b (r/R, 



a^ = -—^^p+ ' ' , (10) 

d m r n 

noting that — >■ p near the origin. In central regions a 
modified radial coordinate, 2 = r^"'', is useful to soften 
terms involving the density, 

d _ rP d 
dz 1 — p dr ' 



(11) 



The stellar velocity dispersions may in principle differ in 
the three orthogonal directions: radial (err), polar (ag) and 
azimuthal (cr^). We will consider spherical models in which 
the transverse velocity dispersions are equal, ao = = a±. 
The radial and transverse velocity dispersions are related by 
defining an anisotropy profile. 



/? = /3(r) 



(12) 



(e.g. Osipkov 1979; Merritt 1985; Mamon & Lokas 2005). 
The radial "pressure" of stars is 



n — 2 



(13) 



and its profile is calculated from a spherical Jeans equation 
for momentum, expressed in a frame co-moving with the 
mean radial velocity {v^) of a mass "packet" of stars (Jeans 
1915; Binney & Tremaine 1987), 



dr 
where g 



gp-k 



2 



2^ 



P-k- 



dvi^ 

dt 



(14) 



■'G{m^ -f ms)lr'^ is the gravitational field. 



2.3 acoustic & orbital timescales 

The acoustic crossing time of the dark matter from the origin 
to the outer surface is obtained by radial integration of 

did ^ 1 
dr , Arr? 



(15) 



For the stars, the analogous "acoustic" time is: 



dr 



dpk 
dPk 



dp* 

dr , 



'APk 

dr 



(16) 



Assuming stationarity, then substituting (10) and (14) sim- 
as the factor cancels the singular part of p^ in any ODEs. plifies (16): e.g. if /3 = then dtk/dr — \J ak/gof'o- From 
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dimensional considerations, the crossing time should scale 
with the halo radius as td ~ R*^^. Appendix B shows ex- 
amples of the exact relations. When calculating the natural 
modes of the galaxy, we should expect oscillations to scale 
in proportion to the angular frequencies ^Id ~ 2-K/td and 
0.^, = 27r/i^. Once the eigenfrequencies (a;) of a galaxy are 
known, one can compare them to local keplerian frequencies, 



G(mTir + TTld) 



(17) 



In regions where k ^ uj, perturbations should only affect 
stellar orbits adiabatically. Wherever n <^ u), tidal shocking 
may occur when disturbances exceed the linear regime. 



2.4 perturbative displacements and gravitation 

For each of the variables describing local properties of the 
dark matter or stars, we define a first order perturbation. 
For instance for a density variable we have 

p = po + pee'"' ^ po(l + \p ee''^*^ (18) 

where the subscript "0" denotes the stationary solution; 
the "~" part describes the response to perturbation; the 
A function describes the perturbation in relative terms (e.g. 
Ap = p/po). Perturbed variables such as Ap and p are com- 
plex. The small constant e is the scale of perturbation; 
Higher order terms are discarded (e^ — > 0). The frequency 
of perturbation is complex (lj = uj^ + iij^ ) with the real part 
describing oscillations and the imaginary part describing the 
tendency to growth or decay. 

As in Eddington (1918) we shall consider lagrangian 
perturbations. For a mass shell that sits at radius r,, in 
equilibrium, we write its displacement and velocity under 
perturbation: 



r — rg + r ee 

~ iu}t 

V = tiur ee 



(19) 
(20) 



The dark matter and stars at any r„ may be displaced dif- 
ferently, so we treat them in two slightly different lagrangian 
frames, and distinguish their coordinates carefully: 



rd^ra+ Td ee 



ro I 1 + A,.^ ee 



iujt • \ iuit 

= lojr^ ee = iluTq Ar^ ee 



(21) 
(22) 

(23) 
(24) 



As the perturbed variables are complex, the mass con- 
stituents may be displaced with different phases (i.e. it is 
possible that arg(rd) 7^ arg(ryf)). Time derivatives are 
straightforward, e.g. for density. 



— _ itot ■ iuit 

-pee = iLjpoee 



(25) 



dt dt 

Expressions for the perturbation of each radial gradient term 
depend on whether we consider the comoving frame of the 
stars or the dark matter. Let's abbreviate and "— " for 
and "d" subscripts. Mass conservation requires that 
p±r±dr± = pipr^drp, and then in the ± frame we have the 
operators 



d^ = d^,+^^^^+^^^^^'' d^ 



(26) 



The gravitational field provides the only coupling between 
components. Consistent evaluation of its perturbations (in 
both frames) requires special care. The stellar and dark 
shells that rest at in the stationary solution will feel the 
same field strength there {go) but the perturbed versions of 
g differ since the shells oscillate through displacements with 
different phase and amplitude. Define 



g± = 90 + g± ee 



(27) 



For newtonian gravity, we have go — —Gmo/r^ with mo — 
mo{< r^) = mo_ + mo^ being the total mass contained 
within in equilibrium. Poisson's equation applies at all 
orders of perturbation, which implies 



^ = -^-4.G(p, + ,-) 
dr± r± 



(28) 



Combining (26) and (27) gives the perturbed part of g± as 
felt by the ± component. Denoting g± = —Gm±/r'^ gives 

^dm± ^ -2goro (Ap^ -I- 3Ar±) 



— ATvGrl [2po±Xr± + P0:p (Ap^ — Ap^ + 2\r±)] , 

(29) 

which has the advantage of being nonsingular and m-i- = 
m_ = at the origin. 



2.5 dark matter perturbations 

The dark matter perturbations are considered adiabatic (s 
is spatially and temporally constant). The equation of state 
links the pressure and density perturbations: 



^Pd = -H- = 7 >Pd = 7 — ■ 
Mass conservation then yields 



dr. 



— — Ap,-2A.,^ - — 



P. r 
2- 

r, 



dXr 



1 /A 



dr„ 



3\r 



(30) 

(31) 
(32) 



The lagrangian equation for momentum (3) decomposes un- 
der (26) to give the hydrostatic condition. 



dPod 
dr„ 



'POd rgXr^ - 2gQ\r^+gd) = podSo Ap^ , 



and a perturbation equation, writable in two forms: 

d^ 
dr-Q 

dAPd _ 30 , . 

-^^--^lA^'d - APd) 

'^^Q "Od 

with the abbreviation 

aj^r-r,Arj — 2goArj 



(33) 



(34) 



(35) 



Apd = 



30 



^ = ^- 2 + ^i^)A,,.(36) 
mo V K-^ 



The presence of the stars is felt through the gravitational 
perturbation obtained from (29). 
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2.6 stellar perturbations 



2.7 boundary conditions 



For the collisionless stellar component, the perturbed mass 
conservation equation resembles the dark matter version, 
but without an equation of state to link density and pressure. 



dr^_ 


P-k 










dro 


' 



- 2- 



' 

Ap* + SAr^,) . 



(37) 
(38) 



Next, we define perturbed forms of the radial and trans- 
verse velocity dispersions of the stars. 



2 2,9 iujt 

a,. — ag + cr^ ee 



(Tx = (1 — /3)o"o + ee'' 



(39) 
(40) 



where and are the perturbations of the radial and 
transverse velocity dispersions respectively. The stationary 
models in this paper assume isotropy (/3 = everywhere) 
although orbital anisotropics occur transiently during per- 
turbations. Nonetheless let's retain /3 terms in the general 
formulation. If the perturbations are radial then there is no 
torque. If the specific angular momentum of a mass element 
is conserved during a radial displacement then each local 
value of r^a']_ stays constant, eliminating a variable: 



2(1-/3)7. 2 , .2 



(41) 



Expressing the radial pressure perturbation as P = pao + 
po<y'^, provides a convenient relation, 



(^0 (Apy — Ap 



(42) 



The Jeans equation (14) separates into stationary and 
perturbed parts. 



dPo, 



2/3P,u 



dXp^ ('^^r^+fl* - 2goAr^ - ,goAp^) 



^-[(2-5/3)A.^+(l-/?)Ap^-Ap^] 



(43) 



(44) 



Given another Jeans equation for evolution of the radial 
stellar pressure (Appendix C), 



0: 



dP^ 



, 2P*"* / 2 , r, 2 2 



(45) 



and performing the perturbation expansion yields an equa- 
tion in the stationary variables (which reduces to a tautol- 
ogy, = 0) plus an equation in perturbed variables. The 
latter, first order equation is: 

~ dr . 

= icjP^ -f 3ici;Po^ -f 2ia;(l /3)P(,^ Ar^ . (46) 
drf, 

Substituting (37) into (46) leads to an algebraic relation, 
= Ap^ - 3Ap^ + (2/3 - 4)A,^ (47) 
which eliminates the Ap^ terms in (37), (38) and (44). 



At the origin, by symmetry, there cannot be any oscillatory 
displacement of mass, and therefore 



lim r . = lim r . 







(48) 



In relative terms the offset is free; any small nonzero value 
of Arj is acceptable at the inner boundary. Regularity of the 
mass equations (32) and (38) then gives: 



lim Ap = — 3Ar . ; 



Ap^ = -3jXr^ . 



Combining (47) and (49) implies 
^lim^Ap^ = -(5-f 2/?)A,^ , 



(49) 
(50) 

(51) 



and P^ — Po^^Xp^ using the central pressure from the sta- 
tionary model. Analogously, the central pressure oscillation 
of dark matter is Pd — Po^Xp^. Assuming that the stars 
and dark matter are comoving in the dense inner regions 
requires Ar^ = Xr^ at the origin. This completes the set of 
inner boundary conditions. 

As (TOj — at the halo surface, (35) usually tends to- 
wards singularity. To avoid this, the natural modes must 
satisfy an outer boundary condition: 



lim, (APd 



Apj=0 



(52) 



Condition (52) is only satisfied at uj eigenfrequencies. The 
stellar mass profile lacks any finite outer boundary. There- 
fore there are no outer boundary conditions applicable to 
the perturbed variables. 



2.8 numerical integration and solution 

To calculate the stationary profiles, trial values are set for 
the pseudo-entropy constant s and the central velocity dis- 
persion of dark matter (To = o"^|r=o)- The system of coupled 
differential equations (4)-(5) is integrated numerically as an 
intial value problem in radius (using the rk4imp and rkSpd 
integrators from the Gnu Scientific Library^). In practice it 
is necessary to change the independent variable (the inte- 
gration coordinate) partway through the integral. Near the 
origin a softened radial coordinate (11) improves accuracy 
in the stellar cusp. In the outskirts it is more practical to 
use o"! or ln{a^) as the independent variable, and integrate 
to a standard limit at — 10~* or 10~^^. (The chosen limit 
doesn't change any results, as long as it is consistent across 
all runs.) At the outer boundary (radius R) the total dark 
mass is recorded (M^ = m^\r=R). 

At a higher level of the program, an amoeba search 
routine reiterates this integration for different trial values of 
{s,To), seeking to minimise the difference between the ob- 
tained {R, Md) and target values chosen by the user. Usually 
these targets are set to collect a sequence of solutions of dif- 
ferent R but equal Md, i.e. describing galaxies of identical 
mass composition but different compactness. 



http : //www . gnu . org/sof tware/gsl/ 
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Once the {s,To) values are known for a desired sta- 
tionary model, the stellar pressure profile (43) is integrated 
inwards from infinity and tabulated finely. Trial values are 
chosen for the complex frequency of perturbation, uj. The 
stationary model's ODEs are integrated again simultane- 
ously with the coupled ODEs for the perturbed variables: 
(29), (31), (32), (34), (35), (37), (38), (43) and (44). When 
the halo surface is reached (r — >■ R) the outer boundary con- 
dition (52) is tested there, to assess the closeness of lu to an 
eigenvalue. Using a complex test score. 



Z = Zb. + iZi = (Apj 



(53) 



a root-finding routine refines the uj trial values until Z 0. 
To seek eigenfrequencies on the real and imaginary axes of u, 
a bisection root-finder is sufficient. To find eigenfrequencies 
elsewhere in the complex plane, a two-dimensional amoeba 
searcher is unleashed on the Z = Z{uj) landscape. Due to 
symmetries of the ODEs in their oj^ factors, the eigenfre- 
quencies can be: purely imaginary (§3.1); purely real (§3.2); 
or occuring in complex conjugate pairs (§3.3). 

If numerical integration traverses many orders of mag- 
nitude in r then unfortunately there are instances in which 
one or more of the A variables overflows or underflows the 
computer variables. In practice, numerical integration to- 
wards the outer boundary produces diverging |Apj|, if uj 
is poorly chosen. This must be tackled through a practical 
numerical trick. One can pick the worst offending variable 
(Aj — Aj^+iAj) and assign its magnitude to a new dependent 
variable, ^ = (A^ -I- A^)^''^. The code tracks the logarithmic 
profile of this quantity, by numerically integrating an extra 
ODE simultaneously with the other ODEs: 



dlnC 



A ^i^ + A ^ 
dr A2 -f A2 V " dr Mr 



(54) 



All the A variables can then be expressed in a locally nor- 
malised form, Afe = Afe/^. If each corresponding ODE (fc) is 
abbreviated as a function dAfc/dr = /ifc(Ai, A2, ...) then 



_d_- _ fafc(Ai, A2, ...) _ Afc ding 
dr g g dr 



■ /ifc(Ai, A2, ...) — Afc 



ding 
dr 



(55) 



The simplification hk(X)/^ = /ifc(A) depends on the linearity 
of the perturbations. Evolving (54) for g = g(r) in conjunc- 
tion with the ODEs of the modified form (55) keeps all the 
A variables normalised to moderate values that do not over- 
flow the floating-point representations. Keeping the running 
values small also avoids driving the numerical integrators 
into stiff or unreliable domains. Near the outer boundary, 
where Ap^ tends to become large, the running values of Ap^ 
and Apj can be kept on the order of unity or smaller. This 
prevents numerical noise in the sign of a score Z = Z/g that 
is based on (53). For presentational purposes, the A profiles 
can be rescaled by g afterwards in post-processing. 

For brevity and computational ease, we choose physical 
units such that G = 1, along with the stellar half-light radius 
Rs = 1 and the stellar density there pc — 1. Results can be 
scaled into metric or astronomical units by inserting fiducial 
values of these constants for a known galaxy. 

Families of models will be labelled by their total stellar 
mass fractions, J-V = i\f-*-/(Afj + M-^)- This paper presents 




-0.2 0.0 0.2 

angular frequency u'/\u\a. 



0.4 



Figure 2. Plot of the difference between calculated and ideal 
values of Ap^ at the outer boundary, for trial uj values. The hori- 
zontal axis combines real and imaginary oj adjacently. Sharp dips 
reveal eigenfrequencies. For the black curve, the stars were a fixed 
background potential. For the orange curve, the stars were free 
to oscillate. The "growing mode" dip disappears when stars and 
dark matter are gravitationally coupled and both free to move. 
This galaxy has parameters F = 8, T-t^ ~ 0.163, R = 100 Re- 



galaxy models for a variety of choices of J^^^ . A fiducial "cos- 
mic" baryon fraction ~ 0.163 is defined as in Saxton & Wu 
(2008) and Saxton & Ferreras (2010). For poorly understood 
reasons, observed densities of baryons (e.g. Persic & Salucci 
1992; Bell et al. 2003; Xue et al. 2008; Anderson & Bregman 
2010; McGaugh et al. 2010) and dark matter (Karachent- 
sev 2012) differ from cosmic mean estimates (except in rich 
galaxy clusters). There may be an unseen ambient sea of 
both materials, free from self-bound objects. Matching the 
abundances of observed galaxies to simulated CDM halo 
distributions implies baryon fractions < 0.04 depending on 
galaxy mass (e.g. Guo et al. 2010; Papastergis et al. 2012). 
This may mean that galaxies lost most of their primordial 
gas endowments, or it may indicate that real halo popula- 
tions differ from CDM. Cosmological initial conditions are 
beyond the scope of the present analyses, so J^.^ is for now a 
free parameter. Nonetheless §3.1 shows that spherical SIDM 
galaxies have stability limits that may influence the emer- 
gent distribution of J-".^ . 

In coUisionless CDM, the outer halo density drops grad- 
ually without a well-defined boundary. These simulated ob- 
jects are conventionally circumscribed by the "virial radius," 
containing some multiple of the cosmic critical density, e.g. 
Pv = lOOpcrit ~ 1.4 X lO*m0kpc~^. If the natural density 
and radius scaling of SIDM haloes is not too different from 
CDM then we can estimate upper limits to the halo radii of 
represenative galaxies, R'^ < 3Mi^ /AirpvJ-i^. For the exam- 
ple of M87, R/Rc < 30/J■i''^ for M32, R/R^ < 
for a red nugget, R/Rc ^ 1 '"^ ' ttI/^ 

R/Rc < W/T]/"^ (data 
unclear whether this rule-of-thumb applies to SIDM, which 
warrants new cosmic structure formation calculations. 



< 120/J'4 ; ^nd for the Fornax dSph 
Walker et al. 2009). It is however 
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3 RESULTS 

3.1 stability properties 

Classical arguments based on secular conditions suggested 
that F > 6 polytropes should be unstable to runaway ex- 
pansion or collapse (e.g. Ritter 1878; Emden 1907; Chan- 
drasekhar 1939, pp 51-53). This tendency manifests itself 
(in the present formulation) as a conjugate pair of modes 
with imaginary ui values. These "collapse modes" or "grow- 
ing modes" involve no oscillatory motion, only exponential 
growth or decay. In practice such modes are found by search- 
ing along the imaginary axis of the uj Argand plane. 

Next, a polytropic halo can be tested in a static stellar 
background (acting as a fixed potential, with Ar^ = Ap^ = 
fh-i, = 0). The stability properties can be compared to the 
generic system where both the stars and dark matter are 
enabled to move radially. Figure 2 illustrates modes of a 
particular model with F — 8, comparing the results if the 
stars are fixed or free. With stars held fixed, a collapse mode 
appears. With the stars enabled to move, the collapse mode 
vanishes from the model shown. This stabilisation is typical 
of galaxies where the halo extends farther than R > 10i?o. 

The stabilisation can be explained qualitatively as 
follows. SIDM fluid carries pressure waves, and F > 6 
modes are collapsible because pressure may not rise strongly 
enough to counteract compressive perturbations. Stars how- 
ever tend to persist on their orbits unless tidally disturbed, 
and do not conduct true sound waves. They deafen the prop- 
agation of overdensities, and exert a local damping and sta- 
bilising effect. For many galaxy models this suffices to confer 
global stability. However, as we shall see below, in extreme 
cases the stabilisation is insufficient (§3.1.1) or a new type 
of instability emerges (§3.1.2; §3.1.4). 

The rest of this paper assumes that the dark matter 
and stars are gravitationally coupled and both respond to 
radial perturbations. Exhaustive numerical surveys have lo- 
cated the collapse modes as functions of: the equation of 
state parameter (F), the radial extent (R) of the halo, and 
the global richness of stars (J^*). Figure 3 illustrates the 
|w| collapse frequencies as functions of R and F (colours), 
for models where dark matter dominates (small JTy^, in re- 
spective panels). Figure 4 compares the modes when DM is 
present in the "cosmic" abundance, and a case that is DM- 
poor. Realistic haloes could have radii of tens or hundreds 
of kpc, perhaps a range 10 < R/Rc < 10'' depending on the 
galajcy. Ultra-dense and ultra-diffuse cases are also plotted 
to clarify some theoretically notable modes. Two distinct 
types of modes occur in particular domains. 

3.1.1 slow growing modes 

If R is large (a very diffuse halo) and F > 6 then an instabil- 
ity occurs. This is a reemergence of the instability of classic 
one-component polytropes with _F > 6. In these cases the 
stars are unable to suppress instability because the stellar 
density is too low in the outer halo. This happens over a 
larger domain of R if J^^, is lower (dark matter dominates 
stars). As R is increased, the mode turns on suddenly: for R 
below this threshold, a galaxy is safe; for slightly larger R the 
mode operates on its shortest timescale; for very large R the 
eigenfrequency diminishes (|cj| ~ R~^^^). In a survey of con- 
figurations with 7 ^ F s; 9 in the range 0.0001 ^ J"* < 0.64, 



the frequency scale of these "slow" growth modes is lower 
than the effective stellar frequency scale, \uj\ < 0.02ajc. In 
other words, the e-folding duration of the growth / col- 
lapse would last several tens of orbits for most stars. The 
slow mode straddles the accoustic frequency scale: always 
< 8ild for F = 9 haloes of the largest radii considered, 
but < lOd more usually. The |a;|/r2d ratios are lower for 
F = 8 than F = 9, and lowest for F = 7. 

The perturbed mass flux density through a sphere of 
given radius is = iujpof, so the quantity por is a useful 
expression of the radial flows of stars and dark matter. Fig- 
ure 5 depicts the radial flow structure during "slow mode" 
evolution. Amplitudes for stars tend to be highest in the 
centre and decline radially. The dark flow has its minimum 
in the centre. The stars tend to move in the same direc- 
tion everywhere. Dark matter has alternating ingoing and 
outgoing layers, separated by nodes at rest. In the interior, 
dark matter moves the same way as the stars. Then there 
is a node at some radius around > O.lRc to ~ 10i?c (far- 
ther out for larger F). Outside this node, the dark flow is in 
antiphase with the stars. Some models have further nodes 
and flow reversals at larger radii in the outer halo (e.g. two 
far outer nodes in the F = 9 example in Figure 5). The 
stability analyses are symmetric under time reversal, so it is 
impossible to say whether a particular sign implies collapse 
or expulsion of matter. The slow mode could imply contrac- 
tion of all stars and the inner dark matter while an external 
dark layer escapes. It could also mean radial expansion of 
stars and inner halo, while dark outskirts converge inwards. 

3.1.2 fast growing modes 

If R is sufficiently small then one or two growing modes ap- 
pear. They occur for all F choices. The affected models have 
hot and dense haloes, with much of the dark mass residing 
inside the half-light radius. These "fast" modes occur at high 
frequencies: usually |a;| » uje, (and often \u]\ 3> 10^ for 
larger T^, cases). Where there are two "fast" eigenvalues, 
the uppermost is quick compared to the acoustic crossing of 
the halo {\uj\ > Sid, in many models by a few dex). The lower 
eigenvalue is faster than Sid if F ^ 7, but can be slower than 
fid for F < 7. 

Figure 6 depicts the mass flows due to the lower-cj fast- 
mode for illustrative models with halo radius R = 10-Rc and 
DM enriched over "cosmic" composition. In some ways these 
proflles resemble the waveforms of slow modes: the stellar 
amplitude is high in the interior; the dark matter amplitude 
is peaked farther out. The dark matter has one node near the 
outer boundary, well outside the half-light radius of stars. 
Except in the outskirts, the dark mass is flowing in the same 
direction as the stars. For the corresponding higher-|tj| fast 
modes (Figure 7) there is again just one node in the dark 
matter, but both constituents show a dramatic rise in am- 
plitude towards the outer boundary. The rise is steeper for 
larger \uj\ and larger F. These fast modes act mainly in the 
diffuse outskirts: rapidly driving the stars and dark matter 
in opposite directions. The core is essentially unaffected. 

Since the fast instability occurs for all F, its main causes 
must be traits of the stellar profile. We might suspect the 
stellar cusp and its coldness (dcr^^ /dr > at the origin) how- 
ever tests show fast modes afflicting non-cuspy models with 
a Plummer (1911) proflle of p.^. If the exponential stellar 
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halo outer radius R / 

Figure 3. Frequencies of tiic halo collapse modes as a function of the outer radius of the halo, for models dominated by dark mass. 
Curves are coloured in spectral order to show the SIDM thermal degrees of freedom: from violet {F = 2) to red (_F = 9). Panels portray 
different values of the global stellar mass fraction. As J^^f is increased, "fast" modes rise in frequency, "slow" modes (bottom-right curves) 
recede to affect only the most radially extended haloes. 



fringe is a factor, the fast modes might diflFer in tidally trun- 
cated models (generalising King 1966) which will be anal- 
ysed in a future paper. In the present model, the index of the 
stellar gravitational field (d In gi^ /d In r) curves significantly 
near r ~ 107?e but settles to the asymptote —2 farther out. 
If the halo surface occurs at r < 30-Rc then the oscillatory 
restoring force varies in form during the cycle. Fast instabil- 
ities may result from the inadequacy of restoring forces for 
haloes that truncate in this region. 



3.1.3 safe configurations 

Fast modes and slow modes are absent throughout a special 
intermediate range of R. In these cases the interaction be- 
tween stars and the halo suppresses the collapse mode (e.g. 
in the "coupled" case of Figure 2). These haloes with > 6 
are stabilised by the mingling of coUisionless stars, contrary 
to naive textbook expectations. The three radial domains 
— defined by the "fast" instability, stable zone, and "slow" 
instability — imply limits on the compactness of the sta- 
ble endpoints of galaxy evolution. The fast growing modes 
mean that models are unstable unless R exceeds some lower 
limit (which is a slowly varing function of F and J-V)- The 
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halo outer radius R / 

Figure 4. Frequencies of tiie halo collapse modes as a function 
halo radius R for galaxies with cosmic abundances (top) and DM- 
depleted cases (bottom). The halo F values are coloured as in 
Figure 3. As is increased, the "fast" instability mode begins 
to affect radially large _F > 6 models as well as the compact and 
-F ^ 6 models. These modes become ubiquitous for all ij, if F > 6 
and stellar mass is sufficiently dominant. 



slow mode only afflicts F > 6 cases, and it applies an upper 
limit on R. This upper limit is a sensitive function of J^^. 
For models with F > Q that are poor in stars {T-t, < 0.001), 
the slow- and fast-mode domains leave only a narrow band 
of safe R. This safe band becomes wider as increases. A 
numerical survey indicates that the slow mode (upper limit) 
vanishes for sufficiently high F^. (It persists till T^, < 0.0004 
for F = 7, J"* < 0.004 for = 8, and J"* ^ 0.01 for F = 9.) 
From this point till J^-^ lS 0-4, galaxies with arbitrarily wide 
outer radii R are stable. Figure 8 illustrates the stability lim- 
its on the central dark density (p;j|r=o) as a function of J^^. 
Figure 9 shows the corresponding limits on the dark mass 
fraction within the half-light radius (/d(< l^c))- The latter 
quantity can be measured through gravitational lensing and 
measurements of kinematic tracers. Significantly, it appears 
that a SIDM halo around an isolated de Vaucouleurs type 
galaxy is only stable if the dark mass within IRc is less than 
« 90% of the total, or less depending on J^^ and F. 

Mass modelling of observed elliptical galaxies often in- 
dicates /d(< iRc) values in the range of a few tens of per- 
cent (Loewenstein & White 1999; Ferreras, Saha & Williams 
2005; Thomas et al. 2005, 2007; CappeUari et al. 2006; 
Bohon et al. 2008; Tortora et al. 2009; Weijmans et al. 2009; 
Saxton & Ferreras 2010; Grillo 2010; Memola, Salucci & 
Babic 2011; Murphy, Gebhardt & Adams 2011; Norris et al. 
2012; Tortora et al. 2012; Wegner et al. 2012). If these galax- 
ies are stable within polytropic haloes, we could broadly in- 
fer 0.005 < J-^f < 0.3 independently of F. Observational 
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Figure 5. Radial profiles of the slow mode in DM-rich galaxies 
(T^ = 0.0001) with large halo (R = 9000/?o) and F = 7,8,9. 
The horizontal axis is radius (r); the vertical axis is a density- 
weighted displacement or flow, logj^Q(pr). The stellar profile is 
black ("-i^-"); dark matter is coloured ("d", using the palette in 
Figure 3). Curves are solid or dot- faded in the parts where the sign 
is positive or negative respectively. Downward spikes are nodes at 
rest. The vertical normalisation is arbitrary. 

estimates of the stellar mass are consistent with this range. 
No large nearby galaxy is yet known to challenge the den- 
sity limit of the fast instability. A few exceptional galaxies 
appear DM-poor at large radii (Mendez et al. 2001; Ro- 
manowsky et al. 2003; Salinas et al. 2012). If their impover- 
ishment is not an illusion due to orbital anisotropics, then 
we might explain them in terms of runaway halo loss via 
either of the instabilities (§4.3). 

3.1.4 fast modes in DM-poor systems 

Entering the domain where stellar mass is globally dominant 
(Figure 4) a new fast mode appears. Beyond some threshold 
(J-^ > 0.4) an island of the fast modes begins to affect large- 
R haloes with F > 6. This implies that DM-poor systems 
have a maximum stable radius, just like the DM-rich systems 
that are subject to slow modes. The u values of these fast 
modes are slowly varying functions of R (horizontal curves 
in bottom panel, Figure 4) which implies they arise from the 
stars rather than the DM crossing time. 

As is increased further, the fast mode affects F > 6 
haloes with all possible R values, and then there are no sta- 
ble configurations with F > 6. For example when = 0.64 
there are no stable models with F = 8 or F = 9. The fast 
modes that afflict these DM-poor models are found to have 




Figure 6. Radial wave profiles of the lower-cj fast growing mode 
in galaxy models with stellar mass fraction = 0.01, halo radius 
R = 10i?o and F = 5,6,8. The horizontal axis is radius (r); the 
vertical axis is logj^g (pr) . Dark matter is coloured ( "d" ) ; the stellar 
profile is black ("■A'")- Curves are solid or fade-dotted in the parts 
where the sign is positive or negative respectively. 

steeply rising waveforms (mainly acting at the halo surface) 
resembling those in Figure 7. Results were not sought be- 
yond J^^ — 0.64, because the system of equations becomes 
stiff in ways that make the uj root-finding slow and arduous. 

3.2 neutrally stable oscillations 

As in single-fluid asteroseismology, (e.g. Eddington 1918; 
Cox 1980; Pijpers 2006; Aerts, Christensen-Dalsgaard & 
Kurtz 2010) each polytropic halo exhibits an infinite se- 
quence of radial pulsation modes at discrete real values of 
oj. These radial modes are neutrally stable: neither growing 
nor decaying. This is to be expected, as the present analysis 
is purely adiabatic; omitting any dissipation or mixing pro- 
cesses (e.g. shocks in the dark matter; violent relax;ation of 
the stars; chaotic migration of stars). 

3. 2. 1 flow profiles of standing waves 

Figure 10 illustrates examples of eigenfunction waveforms 
expressed in terms of oscillating mass flux density, (pf ) for 
the dark matter (coloured curves) and stars (black curves). 
For the stars, flow amplitudes are higher in central regions 
and decrease with radius. For the dark matter, the peak am- 
plitudes occur at intermediate radii. Such profiles are generic 
to the radial pulsations of all models with the [F, F^, R) val- 
ues inspected so far. Near the centre, the dark matter and 



Figure 7. Wave profiles of the higher-w fast growing mode cor- 
responding to Figure 6. The displacement amplitudes rise steeply 
with radius; amplitudes are most significant at the surface and 
negligible in central regions. For larger F the ampitude rise is 
greater and the DM node is at smaller radius. 

stars oscillate in phase, but in outer regions they can oscil- 
late in anti-phase. 

For each constituent there are nodes at specific radii, 
where the matter remains locally at rest, and where there is 
a local reversal in the flow direction. The number of nodes of 
the dark matter oscillations is a harmonic number, A'^ ^ 0. 
Values of N uniquely label the eigenmodes. As N incre- 
ments, additional nodes appear in the outer fringe of the halo 
(e.g. between 50i?G and 100i?o for the A = 1, ... 4 modes in 
Figure 10). Most of the nodes appear far enough out that 
the keplerian orbital frequency is locally less than the pul- 
sation frequency: k < (marked with arrows in Figure 10). 
If wave amplitudes grew to the nonlinear regime, some tidal 
shocking and orbital migration of stars could occur. Such 
dissipative effects scale as oc /9V • g which is of non-linear 
order ~ , while the radial displacements are of order ~ e. 
Inwards/outwards pulsational displacements are an immedi- 
ate effect, but orbital heating (and net outward expansion) 
could become significant gradually as more cycles elapse. 

Depending on A and the galaxy model, a number of 
nodes can also occur in the stellar profile (A.^). For a 
galaxy with F = 9, J-V = 0.163 and R = 100-Ro, the node 
numbers are (A, A*) = (1, 0), (2, 0), (3, 1), (4, 1), (5, 1), (6, 2) 
. . . (9, 3), . . . (18, 5) and so on. Figure 11 shows the dark and 
stellar node positions for the first thirty oscillatory modes 
of this model, and for equivalent models where F = 3 and 
F = 6. It appears that A^ never decreases as A increases, 
but there is no explicit link between these numbers. The se- 
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Figure 8. Central density of dark matter at the limits of fast 
and slow instabilities (solid and dashed lines respectively). The 
domain between the curves is safe. The horizontal axis is the 
global stellar mass fraction 



quence of (A'', A*'.^ ) and their radii provide a fingerprint of a 
particular galaxy, and must be obtained numerically. 

Figure 11 also shows that the dark nodes are more 
crowded near the outer boundary (i?) in halo models with 
larger F. This is a consequence of slower wave propagation 
in the fringes, since small-cored (large F) models have pro- 
portionally lower G profiles in the outskirts. Stellar nodes 
are more sublty affected by F. For models with small F, 
the stellar nodes are more numerous and reach smaller 
radii for given A''. Haloes with large F tend to have fewer 
stellar nodes, and they are distributed farther out. Vary- 
ing the stellar mass fraction has slight effects on the posi- 
tions of the nodes in the outskirts. Across a wide domain 
0.002 ^ s; 0.48 (but with fixed N, F and R), the stel- 
lar nodes do not shift significantly. The outer nodes of the 
dark matter appear more affected: they occur at radii a few 
percent farther out if is very small. 

For waves that have at least one DM node (A' ^ 1), 
the radius of the innermost node varies more widely than 
any other node: from ~ 0.1i?c to ~ 10-Rc in the domain of 
0.002 s; J"* s; 0.48, 2 ^ 9 and 10i?c < -R 100i?c. Its 
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Figure 9. Stability limits on the dark mass fraction within li?o 
corresponding to states in Figure 8. The grey dot-dashed line 
indicates the maximum dark mass fraction possible for given T^, . 
Estimates for individual elliptical galaxies are black dots with JT^ 
upper limits (Thomas et al. 2005; Murphy, Gcbhardt & Adams 
2011; Norris et al. 2012). Thick arrows show /^(l-Re) ranges from 
galaxy samples (Bolton et al. 2008; Thomas et al. 2007; Tortora 
et al. 2012) drawn here at arbitrary T.),. 



location varies only sublty with harmonic number N . Fig- 
ure 12 depicts the detailed distribution of nodes in the in- 
ner, astronomically observable regions (r < 20i?e) of galax- 
ies with _F = 3 and _F = 9 where the halo radius is 
R/Rc ~ 10, 100 and the DM content is rich or "cosmic" 
(7v^ — 0.002, 0.163). All else being equal, the innermost dark 
node tends to occur at smaller radii if the galaxy is DM-poor 
(large and farther out if the galaxy is DM-rich (small 

J-ir). For DM-rich haloes, a larger radius R shifts the in- 
ner dark node to larger radii. For some DM-poorer haloes, 
larger R shifts the inner dark node inwards. Comparison of 
the F = 3,9 cases in Figure 12 also shows how the lower- _F 
models have more stellar nodes at radii small enough that 
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Figure 11. Node locations of dark matter (• coloured) and stars (i^ symbols) for the TV ^ 30 radial pulsation modes of galaxy models 
with = 0.163, R = 100i?c and _F = 3, 6, 9 (left, middle and right panels respectively). The black curve marks the radius where = lj. 
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Figure 12. Node locations for Af < 30 pulsations in the inner regions of model galaxies with _F = 3 (top row) and F = 9 (bottom row) 
and various halo radii and stellar mass fractions (_R = 10, 100-Rc and J^^ = 0.002,0.163 as annotated). The innermost dark node occurs 
at a smaller radius if the stellar mass is relatively large (big J^(^). Symbols and the k = uj line arc marked as in Figure 11. 



ripples might be observable in the profiles of starlight and 
stellar velocities. 



3.2.2 frequency spectrum 

Figure 13 illustrates the oscillatory spectra for various 
baryon fractions (J^^^.) and thermal degrees of freedom F = 3 
or F = 9. Among these and other results, it is usually 
found that most eigenfrequencies cj are approximately lin- 
early spaced in relation to N/t^ (where is the acoustic 
crossing time). Typically Ai/ ~ 0.5 /td (or Auj ~ 0.5 fid) 
independently of (F,Fi,,R). If the mode spacing of a real 



galaxy halo were measurable, this information would pro- 
vide a diagnostic of the halo radius via its crossing time (see 
Appendix B). The confirmed detection or inference of the 
lowest mode would constrain (F, F^, R) more tightly. 

There are domains where the lowest mode deviates from 
the linear spectrum. For galaxies with sufficiently compact 
halo (small R) there exists a fundamental oscillation that 
has no dark matter nodes in the waveform (A'' = 0). Keep- 
ing the same {F,J^ir) but increasing the halo radius (71), 
there is a narrow range of R where the fundamental fre- 
quency detaches from the rest of the spectrum, and drops 
to zero frequency (loq — 0). For radially larger haloes, the 
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Figure 10. Radial wave profiles of tlie first four oscillatory modes 
in a galaxy model with _F = 9 halo of radius 100-Rc, and JT^ 
0.163. The vertical axis is logjQ(pF). Dark matter is coloured; the 
stellar profile is black. Curves are solid or faded in the parts where 
the sign is positive or negative respectively. Arrows indicate the 
radius beyond which the keplerian frequency is less than the radial 
oscillation frequency (k < oj). 



lowest mode is the first overtone (A'^ = 1). For F < 6 the 
fundamental mode drops out somewhere from a few 0.1i?o 
to a few 10 -Ro. The dropout threshold R is greater in DM- 
rich systems (small J-'*). For haloes with F ^ 6 the A'' = 
dropout occurs below R < O.lRc, so that all such galaxies 
of realistic size lack a fundamental mode (if F ^ 6). In the 
(F, JV) = (9,0.002) and (9,0.010) cases of Figure 13, we see 
a further dropout of the N = 1 mode when R is very large. 
The dropout range in R is not in any obvious way related to 
the onset of fast and slow instabilities. Apart from influenc- 
ing the dropout radius, the stellar mass fraction (J-V) hs'S 
little effect on the shape of the oscillatory spectrum. 

For halo sizes that are realistic for elliptical galaxies 
(10^ < R/Re < 10'^) the fundamental mode is nonexistant. 
If so then all oscillatory modes have at least one dark node. 
The innermost dark node turns out to occur at radii suffi- 
ciently small that there may be consequences for the visible 



galaxy. Since this node location is nearly independent of A'^, 
it seems possible that dark eigenmodes could excite each 
other via disturbances in the DM close to this region. 

The stellar nodes are stationary points in the radial 
displacement of the stars. The mean radial velocities would 
be zero at these points of rest, and opposite signs on either 
side of each node. Rippling overdensities and underdensities 
of stars would appear between consecutive nodes. Such den- 
sity variations might be observable if any stellar nodes occur 
within r < lOiic, where the surface brightness is practically 
measurable. In the present model calculations, stellar nodes 
tend to occur farther out than the dark nodes, but some stel- 
lar nodes can appear at smaller radii if F is small or if the 
oscillations are of high order. For observable ripples to ap- 
pear, the order of the overtone needs to be large (A^ ^ 30). 

3.3 mixed complex modes 

As emphasised in §2.8, the A ODEs only depend on the 
perturbation frequency via oj'^ terms. For a one-component 
polytrope the eigenfrequencies should be purely real or 
purely imaginary. The two-component galaxy model how- 
ever is an innately more complicated (double) oscillator, due 
to the g± gravitational couplings. Even if Arj and Xr^ are 
in phase at the inner boundary, they may be out of phase 
at other radii due to the different respective wave speeds. It 
is possible for some eigenfrequencies to occur off the real 
and imaginary axes, in conjugate pairs (oj = uj^ + ioj^, 
uj* = uj^ — iojj). 

In models where such modes do occur, they coexist with 
the ordinary oscillations and growth modes. Finding them 
requires extra effort to explore the ui plane in 2D. The land- 
scape of outer boundary scores, Z = Ap^ — Ap^, exhibits 
steep gradients [dZ/du), which easily confounds "amoeba" 
or "steepest descent" solvers. It is more practical to map Z 
on a fine rectangular grid. The desired solutions appear as 
contour rings or spikes on surface plots of log \ Z\ (e.g. as in 
analyses of shock instabilities Toth & Draine 1993; Saxton 
et al. 1998; Saxton & Wu 2001) or by tracing the intersec- 
tions of Z-^ — and Z^ — contours. When an eigenfre- 
quency is approximately located, an amoeba solver can be 
unleashed in that vicinity to refine the ui solution. 

Figure 14 locates the first few complex modes of a rep- 
resentative system. Complex modes occur in a minority of 
the galaxy models investigated so far. These tend to have 
radially large enough {R > 10^ Rc) that the crossing times 
are similar for the dark matter and stars, ti, « td (see curves 
in Figure Bl in the appendices). Such huge halo radii (Mpc 
scale) are astronomically unrealistic, and so the mathemat- 
ical existence of complex modes may be inconsequential to 
galaxies. Complex modes may be of more interest if they 
manifest in other analogous astrophysical systems: e.g. a 
star cluster possessing a gaseous atmosphere; or collision- 
less particles bound inside a fluid star. 

Like the fast growth modes, the complex modes have 
wave amplitudes that are small in the interior and rise faster 
than a power-law towards the halo's surface. Figure 15 shows 
the profile of the mass displacement (pf ) corresponding to 
the lowest-|cj| complex mode of a model with F = 4. As a; is 
complex the eigenfunction has real and imaginary parts too. 
Radial nodes are numerous for all complex modes, because 
^ ild. The innermost node for the dark matter motions 
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Figure 13. Oscillatory spectra for galaxies for various choices of the halo radius R (horizontal axes). Frequencies (vertical axis) are 
scaled by the acoustic crossing time of the halo. Panels show cases where the halo has F = 3 (top row) and F = 9 (bottom row). The 
stellar mass fraction is = 0.002,0.01,0.04,0.163 in the columns from left to right respectively. The frequencies are nearly linearly 
spaced, except in intervals of R where one of the lower modes vanishes. 



is around r ~ 0.13i?o, similar to the innermost node of this 
galaxy's real-oj pulsations. The only complex modes found 
so far are in F = 2,3,4 haloes, where the relatively high 
DM densities in the outskirts enable significant gravitational 
coupling between the components there. For fixed ( , R) 
larger F causes the complex uj eigenfrequencies to be spaced 
farther apart. 

Whether the high-frequency oscillation {ui-^) combined 
with growth (ojj ) is a fatal instability depends on non-linear 
effects. High-frequency rippling ought to induce violent re- 
laxation and orbital mixing of stars. This would introduce 
non-linear damping, and may prevent the mode from grow- 
ing to destructive amplitudes. Proving such damping effects 
would require higher-order perturbative analyses, including 
work terms, (oc pv ■ g ~ e^). 



4 DISCUSSION & IMPLICATIONS 
4.1 reversibility 

The collapse modes of the cores of giant stars are consid- 
ered to cause some classes of supernovae. Are galaxy haloes 
susceptible to dark explosions if they possess the analogous 
collapse mode? The analogy should not be overstretched. 
Stellar core collapse is irreversible due to the accompanying 
nuclear reactions (and loss of neutrinos), and subsequently 
due to an emerging event horizon (e.g. Goldreich & Weber 
1980) . If these microphysical processes have counterparts in 



the SIDM domain then the convulsions of a galaxy-scale 
halo could sometimes have irreversible results. Otherwise, 
the dark matter collapse modes could produce repeated, re- 
versible, non-linear bouncing. At least four phenomena could 
incur some irreversibility: 

(i) If the dark particles can self-scatter, then their self- 
annihilation and decay channels may also be non-negligible 
in some astrophysically attainable conditions (e.g. Buckley 
& Fox 2010; Feng, Kaplinghat & Yu 2010). High densities of 
SIDM might react to yield neutrinos, photons or other es- 
caping byproducts. This is an "unknown unknown," beyond 
the scope of the present paper. 

(ii) Shock-heating of the SIDM would raise the halo en- 
tropy (unless the dark shocks are reversible, as in a scalar 
field model of Peebles 2000) . A shocked halo will eventually 
settle with a different radius R but with [F, Ti,) unchanged. 

(iii) Tidal shocking and mixing of the stellar orbits could 
alter the orbit distribution irreversibly. This could in princi- 
ple affect the evolution of luminosity verses radius relations 
of elliptical galaxies. Early-type galaxies in violently oscil- 
lating haloes might evolve lower stellar densities via tidal 
exchanges driving stellar diffusion radially. 

(iv) Strong compressions of the inner halo could perhaps 
form a central black hole, or feed it significantly. The pos- 
sibility of forming supermassive black holes (SMBH) via 
parsec-scale Jeans-unstable hydrodynamic "dark gulping" 
was discussed by Saxton & Wu (2008). The masses involved 
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are only an inner portion of a halo, and any collapse requires 
-F > 6. Non-linear effects and long free-fall times would pre- 
vent the halo from vanishing into the horizon entirely. 

Ostriker (2000), Hennawi & Ostriker (2002) and Bal- 
berg & Shapiro (2002) also considered growing black holes 
from SIDM, by a mechanism that depends on a gravother- 
mal catastrophe in thermally conductive, weakly scattering 
haloes (not magnetic SIDM, nor scalar field SIDM). Bal- 
berg & Shapiro (2002) predict plausible submassive black 
hole seeds in galaxy haloes, starting from a cored initial 
condition. Ostriker (2000) and Hennawi & Ostriker (2002) 
set initial conditions with a thermal inversion in a central 
density cusp (similar to the hypothetical cusps of collision- 
less CDM). From this they claimed tight upper limits on 
SIDM interactivity. Either way, gravothermal generation of 
a SMBH requires the SIDM scattering crossection to be in 
a special range that ensures kpc-scale mean-free-paths. If 
this requirement is unmet, collapse need not be imminent 
within a Hubble time (Balberg, Shapiro & Inagaki 2002; 
Ahn & Shapiro 2005; Koda & Shapiro 2011). 



4.2 contrast to collisionless models 

To calculate the discrete and continuous modes of entirely 
collisionless galaxies is a complicated task, due to the higher 
dimensionality and freedom available to phase-space distri- 
bution functions; or due to a lack of closure of the Jeans 
equations for large-amplitude perturbations. Elaborate ma- 
trix formulations are normally used (Kalnajs 1977; Poly- 
achenko & Shukhman 1981; Palmer & Papaloizou 1987; 
Weinberg 1991; Saha 1991). These involve choosing orthogo- 
nal basis functions, and practical applications require infinte 
matrices to be truncated approximately. More elegant meth- 
ods exist for special systems (e.g. Sobouti 1984, 1985, 1986). 

Despite the complication of a live halo as a second mass 
component, the present formulation is tractable because of 
several fortunate conditions. The existence of the dark halo's 
boundary conditions constrains the freedom of the collision- 
less stellar component indirectly via their gravitational cou- 
pling and the Ap target. The finite halo radius ensures that 
the collective modes occur at discrete values of uj, rather 
than on continua. The assumption of small purely radial per- 
turbations enables the zero-torque condition (41). (Analysis 
of non-radial modes would require a different approach.) 

Much of the literature concerning the modes of colli- 
sionless galaxies has concentrated on the presence of "radial 
orbit instabilities" which are non-spherical {I — 2) unsta- 
ble modes of systems where stellar orbits are radially biased 
(Polyachenko & Shukhman 1981; Palmer & Papaloizou 1987; 
Weinberg 1991; Berlin et al. 1994; Trenti & Berlin 2006). 
In theory, non-spherical deformations could arise sponta- 
neously in some objects that are stable against radial per- 
turbations. Though this paper focusses on radial modes, we 
can at least say that radially unstable systems are likely 
to be unstable to non-radial motions too. Persistent non- 
radial pulsations are a possibility, and they could alter the 
instantaneous ellipticities and asymmetries of galaxies and 
their haloes. It is at this stage unknowable whether non- 
radial instabilities could destroy radially stable galaxies. A 
general non-radial analysis would be worthwhile, if another 



closure condition could be derived to replace the zero-torque 
assumption. 

4.3 pulsating galaxies 

It was long ago recognised that collisionless galaxy models, 
if built from superpositions of integrable orbits, are capable 
of persistent pulsations (e.g. Louis & Gerhard 1988; Van- 
dervoort 2003). These require sharp gaps in the phase-space 
distributions of stars. Unfortunately the gaps are difficult to 
resolve numerically in A''-body simulations, due to the innate 
limits of mass discretisation. Perhaps because of this sys- 
tematic effect, the possibility of pulsating galaxies is rarely 
mentioned in simulations and observational research. Or- 
bital chaos and diffusion processes might eventually dampen 
these oscillations in real galaxies and star clusters. 

In the present paper, it is shown that galaxies that com- 
bine stars with a fluid-like dark halo can sustain coupled os- 
cillations, without fine requirements on the stellar orbit dis- 
tributions. A disturbed mass of SIDM would echo and jiggle 
with its own internal sound waves. These darkquakes dis- 
place the visible swarm of stars in their orbits. Reciprocally, 
the stars affect the dark waves, e.g. removing the instability 
of many F > 6 models. Concentric radial overdensities of 
apparently displaced stars, if observed in nature, could be 
a signature of galaxy pulsations. If so then these would be 
an interesting diagnostic of the dark physics, enabling stud- 
ies of skotoseismology — a seismology of dark matter. The 
density ripples of high-overtone darkquakes (with kilopar- 
sec wavelengths) would in principle have subtle effects on 
the gravitational lensing profiles of galaxies. 

"Shell galaxies" are observed to have sharp stellar over- 
densities in bands encircling the main host galaxy (e.g. Ma- 
lin & Carter 1980; Quinn 1984; Hernquist & Quinn 1987; 
Prieur 1988; TurnbuU, Bridges & Carter 1999; Wilkinson 
et al. 2000). These are conventionally thought to be phase- 
wrapped remnants of dwarf galaxies destroyed in minor 
mergers. Stars from the disrupted secondary galaxy are 
strewn preferentially in certain orbits about the primary, 
and migrate till they occupy shells (rather than narrow 
streams). This is an attractive model when the observed 
shells differ from the stellar populations of the host galaxy, 
and when shell radii interleave in steps between the left and 
right sides. If however any case were found in which shell 
stars match the host populations exactly, or if any shell en- 
circled the galaxy completely, then we might alternatively 
propose that host stars were displaced radially by dark halo 
pulsations (of high enough order to reach radii < 10-Rc). If 
the shell intervals appear smaller at larger radii, this would 
signify proximity to the halo's truncation edge, R. 

Passive elliptical galaxies observed at high redshifts ap- 
pear to have modern stellar masses (~ 10"mo) but con- 
tained within half-light radii (Re ~ 1 kpc) much smaller 
than modern counterparts (Cimatti et al. 2004; Toft et al. 
2005; Daddi et al. 2005; Trujillo et al. 2006; Longhetti et al. 
2007; van Dokkum et al. 2008; Damjanov et al. 2009; van 
de Sande et al. 2011; Trujillo, Ferreras & de La Rosa 2011; 
Toft et al. 2012; Ferreras et al. 2012). A mysterious process 
drives these "red nuggets" to inflate on Gyr timescales, but 
without appreciable addition of mass, nor any significant 
starbursts. Expulsion of gas by active nuclei or stellar activ- 
ity might conceivably alter the galaxy's potential (e.g. Fan 
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et al. 2008; Damjanov et al. 2009). Gasless minor mergers 
might inflate the fringes of the stellar profile (e.g. Khoch- 
far & Silk 2006; Naab, Johansson & Ostriker 2009). Now 
we might invoke an additional explanation, in terms of halo 
seismology. For stable galaxies, large oscillations could be 
excited by flybys or disturbances beyond the visible stellar 
outskirts. Stochastic excitations might be effective in cosmic 
environments of any density. If higher overtone pulsations 
reach large enough amplitudes, the outer stars may suffer 
tidal shocking or orbital diffusion, leading to the observed 
radial expansion. However if the initial conditions of the red 
nuggets already transgress the stability limits (Figures 8 and 
9) then mass will redistribute spontaneously on timescales 
~ LO^^ , even without external stimulation. Compariing these 
theorised evolutionary channels to empirical redshift-growth 
trends is worth a dedicated study elsewhere. 

4.4 extremes of dark matter content 

The stability limits on spheroidal galaxies with SIDM 
haloes are loosest for stellar mass fractions around J^-^ ~ 
0.04. The constraints are tighter for systems with very small 
or large J-^ (DM-rich and DM-poor extremes). Thus the the- 
ory predicts limits on the stable end-points of galaxy evolu- 
tion. These depend on intrinsic structural properties rather 
than traits passively inherited from cosmological initial con- 
ditions and halo merger histories. 

The strongest finding is the exclusion of DM impover- 
ished galaxies (J^^^. > 0.5) when the halo has a high heat 
capacity {F > 6). If such a system were born, the ubiqui- 
tous fast modes (§3.1.4) would destroy it by one process or 
another. The fiux profiles imply dramatic evolution in the 
outskirts: perhaps runaway loss of the halo (till only a bare 
star cluster remains). Runaway accretion of dark matter is 
another possibility (if an external reservoir is available) until 
J-^^ is lowered to a stable configuration. The observed non- 
existence of stationary isolated galaxy-scale systems with 
■^■k 0.5 is consistent with SIDM predictions. 

At the opposite extreme in DM richness, observed dwarf 
spheroidal galaxies (dSph) are aged, quiescent stellar sys- 
tems with low luminosities 10"^ < L/Lq < 10^ and dark 
masses > 10^ within IRc (e.g. Mateo 1998; Gilmore 
et al. 2007; Strigari et al. 2008; Walker et al. 2009; Walker 
&c Peiiarrubia 2011; Amorisco & Evans 2012; Jardel & Geb- 
hardt 2012). For the Milky Way satellites, 0.80 < /d(< 
l_Re) ^ 0.999 which implies J^y^ < 0.2 for the brightest and 
J^^ < 10~^ for ultra-faint cases. It appears that dSph oc- 
cupy the upper-left, low-J^^^ part of Figures 8 and 9, where 
their stability is problematic. Their stellar density profiles 
are not cuspy like (7) but rather cored like a Plummer (1911) 
model. Preliminary calculations with this profile reveal fast 
and slow instabilities resembling those reported above. Cus- 
piness is unimportant; the effect of tidally truncated out- 
skirts is less clear. It may also be significant that satellite 
galaxies are not isolated (§4.5). Immersion of dSph in the 
host galaxy's halo would reduce the unstable domains. The 
analysis of galaxies with truncated King (1966) profiles and 
a stabilising external pressure awaits a future paper. 

4.5 alternative boundaries and asymmetries 

A real halo need not be isolated in vacuum. External 



medium confinement can stabilise polytropes (e.g. McCrea 
1957; Bonnor 1958; Horedt 1970; Umemura & Ikeuchi 1986). 
This effect may stabilise clustered galaxies and dwarf satel- 
lites that would be unstable in the field. An ambient back- 
ground of unbound dark matter could imply a warm halo 
surface (crj > 0) . If external DM accretes supersonically, sig- 
nals cannot propagate upstream and the isolated halo model 
still holds. If the contact is somehow subsonic then open 
boundary conditions apply instead of (52): waves might es- 
cape rather than reflecting at the surface. If the interface has 
a high contrast (like the solar atmosphere) then the vacuum 
boundary calculations might approximate cj well enough. If 
the halo blends smoothly into its surroundings then the os- 
cillatory frequency intervals might decrease to the external 
Jeans frequency (y/Up^). 

If the outskirts are infinite and fiuid-like but the asymp- 
totic density is near-zero then the eigenmodes would merge 
into a continuum: a functional relation between uj^ and ojj . 
This relation would be found by setting A values at infinity, 
integrating radially inwards and matching the inner bound- 
ary conditions (48)-(51). This top-down approach is used in 
analyses of radiative shocks (e.g. Saxton, Wu & Pongracic 
1997). The stellar profile (7) is radially infinite, though real 
galaxies may trucate at a tidal radius. If so, an outer bound- 
ary condition on Ap^ analogous to (52) would affect the 
modes. If the stellar boundary is cold {a^ — 0) then stel- 
lar nodes would be bunched there, and ripples could appear 
more numerous than in §3.2.1. As long as either the SIDM or 
the stars are spatially finite, a discrete spectrum is expected. 

In weak-SIDM models with effectively coUisionless out- 
skirts (e.g. Yoshida et al. 2000; Vogelsberger, Zavala & Loeb 
2012; Peter et al. 2012; Rocha et al. 2012), the outer halo 
would act as part of the component in this paper. The 
fiuid-like "d" component would reside in isolation within a 
reduced boundary radius R, shortening the crossing time 
and raising the pulsation frequencies. 

Unless a galaxy is well isolated, its halo may receive 
cosmological infall of smaller objects. If the outskirts are 
quasi-coUisionless then these mergers may cause tidal dis- 
turbances that excite oscillations in the core. If the outskirts 
are in a strong-SIDM condition, then weak shocks occur, the 
impactor can ablate and dissolve, and ripples may propagate 
(akin to raindrops in a pond). Mergers and mixing raise the 
total halo entropy (oc Mdlns), which may alter the final 
halo mass-radius trends. 

Non-radial effects could complicate the oscillations and 
stability limits. Disc galaxies (which provide the best evi- 
dence of cores) require axisymmetric models to assess the 
role of rotation (e.g. Hachisu 1986). Tidal fields, cosmic fil- 
aments and mergers could excite non-radial pulsations di- 
rectly. Non-radial instabilities might destroy some configu- 
rations that are radially stable. The directionally dependent 
crossing times of a triaxial halo would split the frequencies, 
giving a non-linear spectrum of Wj modes with angular de- 
pendencies. Since _F > 6 polytropes are less stable than the 
standard F = 3 case, we might expect larger shape distor- 
tions during non-radial pulsation. The greatest attainable 
transient triaxiality must depend on F and the environment; 
the extent is unknowable without 3D simulations. It would 
be interesting to see the results of strong-SIDM F > 6 shape 
investigations extending Peter et al. (2012). 
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4.6 new thermostatistics, dark fluid or exotica? 

Tsallis thermostatistics generalise the Boltzmann entropy 
in an attempt to describe systems where long-range interac- 
tions (such as gravity) are influential (Tsallis 1988; Plastino 
& Plastino 1993; Nunez et al. 2006; Zavala et al. 2006; Vig- 
nat, Plastino & Plastino 2011). If this conjecture applies to 
self-bound bodies of dark matter, then the optimal spheri- 
cal configurations are polytropes, even without short-range 
interactions or scatterings among dark particles. The equi- 
librium solutions for spherical galaxies have the same density 
and potential profiles as for SIDM fiuid polytropes. Saxton 
& Ferreras (2010) fitted polytropic halo models to observed 
elliptical galaxies without distinguishing between SIDM and 
Tsallis scenarios. 

However in their time-dependent and dynamical prop- 
erties, the SIDM polytrope is not entirely equivalent to the 
Tsallis "stellar polytropes" (e.g. Chavanis 2002). Fluids and 
stellar polytropes differ in their phase-space distributions, 
and their microscopic physics. If coUisionless, dark matter 
need not support or transmit pressure waves. For a coUi- 
sionless dark halo (even with a polytropic profile) the nec- 
essary analytic methods would involve phase-space distri- 
bution functions or moments in Jeans equations, with all 
the problems that accompany the treatment of the stars. 
The complexities of two-species coUisionless modelling are 
beyond the scope of the present paper. Such models might 
not exhibit the oscillations and stability properties calcu- 
lated above. Observation of darkquakes would favour SIDM 
over the Tsallis interpretation. 

Besides SIDM, various alternative DM particle theo- 
ries may be seismically distinct. If haloes consist of neutral 
fermions, then degeneracy pressure may explain the cores 
while the fringes are coUisionless (e.g. Munyaneza & Bier- 
mann 2005; Nakajima & Morikawa 2007). Other candidates 
include ultra-light particles with de Broglie lengths reaching 
kpc scales (e.g. Sin 1994; Ji & Sin 1994; Lee & Koh 1996; Hu, 
Barkana & Gruzinov 2000; Rindler-Daller & Shapiro 2012). 
Perturbations and collisions of these "fuzzy DM" galaxies 
produce interference patterns (Gonzalez & Guzman 2011). 
Such exotic halo models deserve their own stability analy- 
ses, in order to predict observable signatures distinguishing 
them from SIDM and CDM. 



5 CONCLUSIONS 

This paper formulates perturbative stability analyses for 
spheroidal self-gravitating systems comprising a fiuid min- 
gled with coUisionless matter. The present application de- 
scribes of spheroidal galaxies comprising a halo of self- 
interating dark matter and an empirical de Vaucouleurs dis- 
tribution of stars (with isotropic orbits in equilibrium) . The 
halo is polytropic, adiabatic, non-singular and radially finite. 
It is isolated, with a vacuum outer boundary. The perturba- 
tions from equilibrium are radial but non-uniform. 

All radially finite halo models exhibit a spectrum of nat- 
ural modes of radial pulsation. These modes are neutrally 
stable, at least in the absence of high-order dissipative ef- 
fects and stellar orbital mixing. The frequencies are almost 
linearly distributed. The frequency interval scales inversely 
with the halo's accoustic crossing time, which depends on 
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Figure 15. Log-amplitudes of density-weighted displacements 
(pf) for stars (black) and dark matter (coloured) in the lowest 
complex mode of the galaxy model with F = A, R = 5000i?,c and 
= 0.5. The eigenfrequency oj = (0.740 -I- 2.751i) X lO'^We. 
Top and bottom panels show the real and imaginary parts. 



the halo radius. The thermal degrees of freedom F inUu- 
ence the frequency gap below the lowest mode. For haloes 
of plausible radius {R ^ 10i?e) and central density, all oscil- 
lations have a dark node within ten stellar half-light radii. If 
detected via its rippling displacement of gas and stars, these 
modes would help diagnose the actual radial extent and to- 
tal mass of the dark halo. Nodes in the stellar oscillation 
are more likely to appear as detectable ripples and velocity 
patterns if the harmonic number A'^ is large, if the halo is 
compact (small R), or if there are few thermal degrees of 
freedom (small F). 

Gravitational coupling of a coUisionless stellar distri- 
bution with a polytropic dark halo significantly affects the 
collective stability of the galaxy. Stability depends on the 
dark heat capacity (via F), the stellar mass fraction (J-*), 
and the halo radius (R). Elliptical galaxies are stable for 
plausible J^^^ and radii (R) smaller than Mpc scale, even in 
cases where the classical pure polytrope would be unstable 
(-F > 6) . This includes the range of dark matter equations of 
state, 7 < -F < 10, favoured by observed galaxy kinematics, 
cluster core sizes and black hole masses (Saxton & Wu 2008; 
Saxton & Ferreras 2010). The dark matter fraction within 
the stellar half-light radius never exceeds /d(< 1-Rc) ^ 0.9 
for any stable de Vaucouleurs type configuration. A spheroid 
born in an unstable state would presumably evolve towards 
one of the stable end-points. 

The stability limits are defined by the onset of three 
types of growing modes. Slow growing modes occur if the 
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Figure 14. Map of the boundary condition score, Z = (Ap^ — Ap), in part of the u) complex plane, for a large- _R galaxy model. Black 
contours show where the real part is zero; violet contours show where the imaginary part is zero. The intersections of black and violet 
contours are complex eigenfrequencies (seven such knots appear in this plot). 



halo is radially large and F > 6. Slow modes resemble non- 
cyclic versions of the pulsation waveforms: maximum stellar 
flows near the centre, and maximum dark flows at interme- 
diate radii, with DM and stars moving in opposite directions 
in some regions. Fast growing modes (quicker than orbital 
periods at f-Rc) afflict galaxies with very centrally concen- 
trated dark matter (regardless of F). fsolated systems that 
are rich in dark matter are only stable if the central density 
is quite low and the halo radius is huge, ff stars dominate the 
galaxy mass, haloes with F > 6 may suffer the fast instabil- 
ity regardless of the radius. Complex modes are a mixture of 
oscillation and exponential growth. They exist in principle 
but they only seem to afflict haloes of unrealistically vast 
radii. Fast growing modes and complex modes entail erup- 
tion or collapse of mass layers mainly in the fringes of the 
halo, leaving the interior less disturbed. 

The identification of these oscillations and instabilities 
(perhaps via visible stellar displacements or ripples in grav- 
itational lens profiles) would place informative constraints 
on fundamental dark matter physics, and distinguish SfDM 
from other candidates, ff proven, more detailed skotoseismic 
analyses would provide complementary probes the mass and 
energy distribution surrounding individual galaxies. This 
model also suggests an alternative interpretive framework 
that may help to explain: the expansion of red nuggets; the 
possible loss of haloes from extreme galaxies; morphologies 
of some shell galaxies; and perhaps some of the drivers of 
early-type galaxy scaling relations. 
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Table Al. Signature radii of perfect polytropic halos with ther- 
mal degrees of freedom F and outer truncation radius R. The 
radii Ri , R2 , R3 and R/^ mark integer steps in the density index. 
The velocity of circular orbits is maximal at Ro- 
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APPENDIX A: CORE SCALES 

A spherical polytropic halo without other gravitating com- 
ponents is a classic Lane-Emden sphere (Lane 1870; Emden 
1907). The core can be defined by locations where the loga- 
rithmic index of density reaches certain values, 

din pd 



Qd 



(Al) 



d In r ' 

with Od = 1 at radius Ri, and Ofd = 2 at _R2. These occur 
in fixed ratios to the halo surface (i?). The index i?3 is im- 
portant because beyond this radius an isolated polytropic 
SIDM halo differs from the NFW coUisionless halo, in which 
Pd ~ at large radii (Navarro, Frenk & White 1996). 
Table Al shows standard radius ratios for various F. The 
radius where the circular orbit velocity is maximal is also 
given in terms of Ro/ R. Other models' key radii were stated 
in Appendix D of Saxton & Wu (2008). 



APPENDIX B: EQUILIBRIUM PROPERTIES 

The acoustic crossing time of the dark matter halo (§2.3) 
scales with the halo radius as td = C* R?^^ , where the factor 
C = C(R,F) varies slowly in R. For a fixed stellar mass 
fraction (J-i^) the factor C is slightly larger for greater F. 
The global mass ratio has a greater effect on C than R or 
F do. The left panels of Figure Bl show the range of val- 
ues obtained for DM-dominated models, and cosmic mean 
composition. The crossing time associated with the stars is 
shorter than for dark matter, for all models except those 
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with very large radii. This is a result of the coldness of the 
dark matter in layers near the halo outer surface, where 
the stellar density profile is assumed to continue outwards 
gradually. Within the half-light radius of the stars, the dark 
matter has a higher velocity dispersion. Consequently, dis- 
turbances traverse the central regions faster via dark matter 
than via stars; and faster via stars in the outskirts. 

The right panels of Figure Bl show the dark mass frac- 
tion within IRc as a function of the outer radius of the halo. 
For models with F < 6, these fractions decline rapidly with 
R. For models with F > 6 the dark mass fraction within 
these central regions is not very sensitive to the outer ex- 
tent of the halo; it is effectively determined by the global 
ratio of dark to stellar matter. In observed systems where 
the dark matter fraction has been estimated from the kine- 
matics of stars, planetary nebulae or globular clusters, the 
range of model /d(< IRc) and /d(< 5Rc) might be used 
graphically to constrain the likeliest {F,J^^,R) parameters 
of the galaxy as a preliminary step to seismic modelling. 



2p 



-{VrVg) - (Vr't'l) + (vf.) + 



{vivo) cot 9 



(C5) 



If we assume that our radial perturbations are small enough 
to preserve a zero skewness of the velocity distribution in re- 
spect to Vr, (e.g. a gaussian form for the stationary model) 
then {v1) — v'^ + a'^ and (Vr) = v{v^ + Scr^). This step is 
similar to the truncation of moment equations in the one- 
dimensional slab models of Louis (1992), but more generally 
allows pressure anisotropy and radial flows. Then (C5) con- 
verts to the "lagrangian frame" (following the local mean 
velocity) by way of 



d d „ 
— = h V ■ V 

dt dt 



to obtain 
dP 



(C6) 



dv 2pv 



° = ^+2^'I + (3P + 2pv^) ^-f ^ (v^ + 2.? - al) .(C7) 



APPENDIX C: JEANS EQUATION FOR 
STELLAR PRESSURE 

The stars are described by a velocity distribution function in 
spherical coordinates, / — f{r,6,cj>,Vr,ve,v^), where mean 
observables at any radius r are obtained by integration over 
velocity-space, {q) = JJJ/gd"^ii. (This appendix will omit 
the -A' subscripts.) In spherical coordinates, the coUisionless 
Boltzmann equation is 

dt ^ dr r 89 r sin 9 d(f> 

+,,|/+,,|/+,,|/ =0 (CI) 

OVr OV0 aVrf, 

with the derivative terms 
r or 

vl cot 9 - VrV0 1 9$ 

^e^— -EF (C3) 

T r o9 

~V4,{vr + vecot9) 1 a* 

T rsui9o4> 

Multiplying (CI) by an observable quantity q and taking 
means over the velocity-space yields moment equations of 
the type found by Jeans (1915). 

Consider the special case where a spherically symmetric 
galaxy is non-rotating: 9 and (f) are irrelevant and the vq and 
velocity components are indistinguishable. Then we have 
f{r,Vr,ve,v^) — f{r,Vr,±ve,±v,j,). If we consider purely 
radial perturbations in which the distribution function re- 
mains even with respect to the transverse velocity compo- 
nents, then some velocity moments vanish, (e.g. (ve) = 0, 
(VrVg) « 0) and {vrvj) = {vrV^) « va\. The gravita- 
tional potential $ is spherical, so d<^/d9 — d^/dcj> = 0, 
and the gravitational field is radial {g = —d^/dr). Under 
these simplifications, multiplying (CI) by and integrat- 
ing over velocity-space obtains a time-dependent spherical 
Jeans equation for the radial pressure, in an eulerian frame: 
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Figure Bl. Scaling properties of solutions for various halo radii {R, horizontal axis), different F values (coloured curves) and stellar 
mass fractions {T^, = 0.002 and 0.163). The left column shows the factor C = t^/R^^'^ in the relation between acoustic crossing time 
and halo outer radius. The black curve shows the crossing timescale of stars (on isotropic orbits) within the halo {t^ replacing t^). The 
right column shows the dark mass fraction within IRc (the half-light radius of the stars), for models of specified halo outer radius R. 



